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KXPLICIT  CHKBYCHEV  ITERATIVE  SOLUTION  OF  NONSELF-ADJOINT 
ELLIPl’IC  EQUATIONS  ON  A VECTOR  COMPUTER 

Introduction 

The  computational  mathematics  literature  is  well  stocked  with 

numerical  methods  for  solution  of  self-adjoint  linear  operator 

equations.  See,  for  example,  the  many  methods  and  references  given 

by  Varga, ^ Birkhoff,^  and  Vichnevetsky However,  many  fewer  approaches 

are  offered  for  nonself  adjoint  problems . One  occasionally  encounters 

the  suggestion  that  the  equation  L(r)^(r)  = S(r)  (with  L the  linear 

operator,  and  S the  known  driving  term)  be  made  self  adjoint  by  an 

t t 

extra  application  of  the  adjoint  of  L:  L = L S . There  may  be  cases 
in  which  this  approach  has  merit.  However,  it  has  two  immediate  draw- 
backs: (H  one  has  to  work  with  a higher  order  equation,  and  (2)  the 

ratio  of  maximum  to  minimum  eigenvalue  amplitude,  an  indicator  of  the 
amount  of  numerical  work  required,  is  squared.  There  are  situations 
in  which  it  may  be  desirable  to  transform  a self-adjoint  problem  into 
a nonself-adjoint  one.  For  example,  in  a heat  conduction  problem,  L 
contains  the  operator  V*KV,  where  K is  the  thermal  conductivity.  If 
K varies  greatly  in  magnitude,  the  eigenvalue  span  of  L will  be  in- 
creased accordingly,  making  the  numerical  solution  more  difficult.  A 
possible  remedy  is  to  divide  through  by  K,  so  that  L contains  the  non- 
self adjoint  operator  V log  K*V. 

An  approach  that  has  gained  recent  attention  is  a combination  of 

approximate  factorization  and  the  conjugate  gradient  method.  See 

Vlchnevetsky,^  p.  60,  or  Meijerink  and  van  der  Vorst."^  With  some 

adaptation,  this  approach  may  hold  promise  for  nonself-adjoint  problems. 
Note;  manuscript  submitted  June  13,  1977. 
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But  the  factorization  stage  is  awkward  for  vector  computers.  One 
must  solve  upper-and  lower-triangular  matrix  equations,  which  resemble 
recursion  formulas.  As  a result,  the  computer  must  initiate  and  com- 
plete certain  calculations  one  at  a time,  rather  than  doing  several  in 
parallel.  Thus  we  wish  to  offer  an  alternative  method  which  allows 
complete  vectorization  (parallel  processing)  for  all  interior  mesh 
points  of  a multidimensional  grid.  This  method  is  a generalization  of 
the  Chebychev  semi -iterative  method  (Varga, ^ Chap.  5)> 

Assumptions 

We  seek  an  iterative  solution  to  the  equation 

U = S,  (1) 

where  S is  a known  source  vector,  and  L is  a matrix  or  finite 
difference  operator.  We  assume  L has  complex  eigenvalues 
with  all  being  of  the  same  sign.  We  also  assume  all  X fall  within 
an  ellipse  in  the  complex  plane  (see  Fig.  l)  whose  major  axis  coincides 
with  the  real  axis.  The  intersections  of  the  ellipse  with  the  real 
axis  are  b - a and  b + a,  with 

b/a  >1.  (2) 

The  semiminor  axis  of  the  ellipse  is  c.  We  require 

c < ( a(  , and  thus  (J ) 

X.  -X.  < X - X . (4) 

i i , r r , 

max  min  max  min 
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We  assume  that  ( 1 ) possesses  an  exact  solution  4 . At  the  end 
of  n iterations,  we  will  have  an  approximate  solution  <^ , whose  error 
is  defined  to  be 

lO  X / \ 

e = (5) 

The  iterative  method  is  to  be  such  that 

e"  = P^(L)£°,  (6) 

where  is  a polynomial  of  degree  n.  Substitution  of  (5)  into  (6) 
gives 

«S"  = P^(L'«5°  - (P^fL';  - 1)$.  (7) 

We  do  not  know  $ in  advance,  but  we  do  know  L 4 = L S for  k > 0. 

Thus  from  (7)  we  must  require  that  the  zero-degree  term  of  Pj^(L)  be  1: 

P„(o)  = 1.  (8) 

We  wish  to  choose  P^  such  that  its  magnitude  is  as  small  as  possible 

everywhere  within  the  ellipse  containing  eigenvalues  of  L. 

Determination  of  P for  a Real  Min-Max  Problem 
n 

The  problem  of  minimizing  the  maximum  value  of  |P^(x)l  subject  to 

P^fO)  = 1 has  a well  known  solution  when  x is  restricted  to  real  values 

between  b - a and  b + a with  b/a  > 1.  The  standard  argument  points  out 

that  the  desired  P is  such  that  all  maxima  of  P have  the  same  value, 
n ' n 

One  immediately  determines  that  P^  is  proportional  to  a Chebychev  polynomial 
The  Chebychev  polynomials  T^  are  such  that 
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and 


T (cos  O' ) = cos  nO', 

n 

(9) 

T (cosh  O' ) = cosh  nff. 
n 

That  cos  na'  is  a polynomial  in  cos  oi  results  from  the  elementary 
identity 


cos(m  + 1 )o  = 2 cos  oi  cos  raO  - cosCm  - l)®.  '' 10  ) 

Equation  i 10  ; is  also  valid  when  cosh  is  substituted  for  cos.  Thus 
if  cos  mo  and  coS'^m  - 1 are  polynomials  in  cos  a,  then  so  is 
cos(m  + l)^.  This  is  the  case  for  m = 0 and  1,  so  it  is  also  true 
for  all  m > 0.  Equation  ( 10 ) also  gives  the  recursion  formula  for 
the  Chebychev  polynomials: 


T_^/x)  = 2 X T (x)  - T /x).  'll) 

mfl  m m-1  ' 

From  (^9)  we  can  see  that  It^  x)l  reaches  a limiting  value  of  1 n + 1 
times  as  x varies  from  -1  to  1.  Thus  the  solution  to  the  min-max 
problem  for  real  x is 


and 

1 P I = ! T 1 ^ ( 13 ) 

n max  n I a y ' ^ ^ 

= 1/cosh/n  cosh  ^ < 1. 
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Determination  of  for  a Complex  Min-Max  Problem 

Slightly  less  well  known  is  the  fact  that  the  Chebychev  poly- 
nomials also  satisfy  the  following  problem:  With 

z = ? + iT)  (14) 

find  such  that  P^  O)  = 1,  and  such  that  every  maximum  of 

jp^iz)!  on  the  ellipse 


has  the  same  value. 

Let  us  now  consider  the  behavior  of  the  Chebychev  polynomials  in 
the  complex  plane.  It  is  convenient  to  express  ? and  Tj  in  terms  of 
dimensionless  variables  , a and  B as  follows: 

T1  = X Tl' 

° (16) 

^ = cos  ot  cosh  S 

1]^  = sin  a sinh  B . 

Here  X^  is  a dimensional  quantity  of  arbitrary  magnitude  having  the 
same  dimensions  as  the  operator  L.  Let  us  first  demonstrate  that 
arbitrary  ? and  TV’can  be  expressed  as  in  ( l6 ) with  real  a and  B. 

One  can  eliminate  cos  a from  ( l6 ) and  obtain  the  following  equation 
for  cosh  B: 
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u(?'2  + n'2  + 1)  + 5'2  = 0, 


(17) 


F(u)  = - 

where  u = cosh^  6.  Incidentally,  one  obtains  the  same  equation  by 
eliminating  cosh  0 and  setting  u = cos^  a.  Note  from  (17)  that 

F(0)  = 

F(  1 ) = - T]  ^ and 
F(?'^  + + 1)  = 

We  see  immediately  that  there  must  be  one  real  root  0 ^ u^,  ^ 1 and 
another  real  root  1^U2^  I + + ~f]'^ . The  lesser  of  the  two 

roots  can  be  taken  to  be  cos^  O',  and  the  greater  to  be  cosh^  6, 
with  a and  0 both  real.  The  signs  of  5 ^ and  T]  can  be  adjusted 
properly  by  proper  choice  of  the  sign  of  0 and  phase  of  a modulo  tt. 
Let  us  return  to  l6)  and  define 

Z = ? + iT] 

= cos  O'  cosh  0 + i sin  <y  sinh  0 (19) 

= cos  {a  + i0  ) • 

From  (9)  and  (19)  we  have 

T^{z' ) = cos(na  + inB  ) (20) 

= cos  00*  cosh  n0  + i sin  no;  sinh  n0 . 

From  (20)  we  find  with  only  a slight  amount  of  algebra 

It^(z^)|^  = cos^  no  + sinh^  n0 . (21 ) 

We  can  see  from  ( l6 ) that  surfaces  of  constant  0 are  ellipses  in 
? ^ - 1|  space : 
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? ^ + (T1  /tanh  B )^  = cosh^  P. 


f22) 


With  S fixed,  oac  can  see  from  ( l6 ) that  the  point  (?',  T] ' ) scans  once 

around  the  ellipse  asa  varies  from  0 to  2tt.  Thus  (21)  shows  that 

)l  reaches  the  same  maximum  value  2 n times  as  z'  scans  around 

the  perimeter  of  the  ellipse  of  constant  B. 

It  remains  only  to  show  that  no  maxima  occur  in  the  Interior  of 

the  ellipse  (22).  We  can  conclude  this  from  the  fact  that  T (z')  has 

n 

n distinct  real  roots: 


Tjz') 


n 

m=l 


- cos 


- 1/2  )tt 


(25) 


where  K is  a normalization  coefficient.  Thus 
n 


T ' z 
n 


)l^ 


n 


cos 


(24) 


This  shows  that  T^  z )]  increases  monotonically  away  from  the  real 
axis  T]  = O),  and  thus  there  are  no  local  maxima.  Thus  all  maxima 
in  the  region  enclosed  by  the  ellipse  (22)  must  occur  on  the  boundary. 
As  pointed  out  by  21  , these  maxima  of  lT^(z')j  all  have  the  same 
value.  Therefore  we  conclude  that  the  Chebychev  polynomials  solve  the 


problem  stated  just  prior  to  (15). 

In  order  to  put  the  ellipse  (15)  into  the  form  (22),  set 


where  f is  a dimensionless  number  yet  to  be  specified.  Substitution 
of  {25)  into  (15^  gives  for  the  perimeter  of  the  ellipse 


5 + (T  ^a/c  )^  = f 


Thus  from  (22)  we  take 


(26) 


B = tanh  ^ c/a 

f = 1/cosh  B (27) 

= (1  - c2/a2)i/2. 

The  desired  polynomial  is  proportional  to  T^(?  + f^l  ),  and  in  fact 

is 


P 

n 


(z  ) 


T 

n 


(28) 


Note  that  in  the  limit  c 0,  we  have  f 1,  and  z | . In  this  limit 
the  complex  solution  (28)  reduces  to  the  real  solution  (12)  as  it 
should . 

What  maximum  amplitude  is  attained  by  28)?  From  (21)  we  find 
for  the  ellipse  of  constant  B 


T (z 
n 


max 


= (1  + sinh^  nB)^/^ 


(29) 


= cosh  nB . 

Substituting  this  into  (28)  gives 
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max 


'50) 


P z 
n 


)! 


T^(cosh  B ) 

T b/a  cosh  B 
n 


where  B is  given  in  (27)-  We  have  dropped  the  minus  sign  from  the 
argument  of  the  denominator  of  (28)  because  T^(xj  is  an  even  or  odd 
function  of  x when  n is  even  or  odd  respectively.  We  know  (50 ) is 
less  than  1 because  the  argument  in  the  denominator  is  larger  than 
that  in  the  numerator  and  T^(x)  is  monotonically  increasing  for  x > 1, 
as  one  can  show  from  (9)- 

Construction  of  the  Iterative  Method 

We  shall  take  the  polynomial  of  (28)  to  be  the  proper  choice  for 
use  in  constructing  the  approximate  solution  to  (1)  according  to  (7^* 
We  must  choose  between  two  approaches  for  generating  the  polynomial 
P (l).  The  first  is  factorization  of  P into  n linear  factors.  We 
have  tried  this  approach  and  found  it  subject  to  roundoff  error 
amplification.  The  reason  for  this  is  that  the  process  is  equivalent 
to  n overrelaxation  steps,  some  of  which  reduce  the  long  wavelength 
errors  at  the  expense  of  increasing  the  short  wavelength  errors.  The 
short  wavelength  errors  are  eventually  brought  back  down,  but  any 
roundoff  error  introduced  while  they  were  large  tends  to  contaminate 
the  final  result. 

The  second  and  more  preferable  approach  to  constructing  Pj^(h) 
is  through  the  use  of  a recursion  formula.  The  only  possible  drawback 
to  this  method  is  that  one  extra  array  of  storage  is  required  for 
retention  of  earlier  iterates.  We  substitute  (28)  into  (7)  to  obtain 
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In  the  numerator  of  this  expression  let  us  express  ^ in  terms  of 
^n  - 1 to  the  recursion  formula  (11),  and  then  express 

Tn  and  ^ in  terms  of  i and  according  to  fjl).  Let  us  detine 

_b 

" fa  (5^) 

= - b/v'a-'  - c^^ 
and  use  = S . 

Then  we  find 

^ T [2T^  q)/faf'(L  - b <5*^  - S i - T^_^(q 


+ $ [1  + (-2  qT^'q)  + V 1 ^ Vl  ^ ’ ^55) 


Notice  that  the  coefficient  of  i is  zero  by  (11).  Recall  that  we 
imposed  this  condition  at  the  outset  in  (8). 

In  order  to  use  this  method  we  must  know  the  first  two  terms  in 
the  recursion.  We  get  these  by  direct  construction  of  P^(l),  and 
obtain 


4’’  = arbitrary  trail  solution  (54) 

4^  = 4'^  - l/b  {h4°  - s) 
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2 T^(q) 

Vi 


((L  - b'^  - S) 


T ,(q)  1 

n-1  /H-l 


/ i-  ^ ^ 

4>  , n > 0 . 


Error  Estimates 


If  B is  an  eigenvector  of  L such  that 


LB  = >B, 


where  X is  a complex  constant,  then 


P (l)b  = P (X)B. 
n'  ' n 


Let  us  imagine  that  the  error  e°  of  (5)  is  expressed  as  a linear  com- 
bination of  the  eigenvectors  of  L.  We  have  assumed  all  eigenvalues 
of  L to  lie  within  an  ellipse  i'15)  in  the  complex  plane.  Then  from 
(6),  (27)  and  (50),  the  coefficients  of  the  eigenvector  expansion  of 
have  each  been  decreased  in  amplitude  by  at  least  a factor 


|Pn(^)l  ^ b^(a,b,c)  = ) /T^  ( b 


In  Fig.  2 we  plot  this  error  limit  as  a function  of  n for  various 


values  of  the  ratio  of  maximum  to  minimum  real  eigenvalue. 
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b+a 

b-a 


'38) 


R = 


} 


and  ellipse  axis  ratio  c/a.  One  notices  that  the  error  limit  decreases 
approximately  exponentially  with  n,  and  that  the  rate  of  decrease 
depends  strongly  on  R.  In  practice,  when  R is  greater  than  about  10®, 
we  use  a filtering  process  to  be  described  later  to  bring  down  the 
effective  value  of  R.  The  slow  convergence  at  large  R is  due  to  the 
fact  that  one  is  attempting  to  reduce  the  error  by  a diffusion  process, 
which  is  slow  for  long  wavelengths. 

Let  us  derive  an  approximate  expression  for  (57)  iri  the  limit  of 
large  R.  It  is  convenient  to  define 


6 = b/a  - 1 

= 2/^R  - 1),  and  (39) 

P = c/a. 

From  (9)  and  the  identity, 

cosh”^  X = log  (x  + - 1 ) {ho) 


we  have 


cosh  (n  log  ) 

E = 

cosh  (n  log  Q^+  n log  Q^) 

= l/[cosh(n  log  Q )+  tanh'n  log  Q ) sinh(n  log  Q_)] 
2 1 ^ 


where 


Q 

1 


1 + P 

Vi  - 


{hi) 


(42) 
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^ _ 1 + 6 + ./2  6 + 6-^  + p-=^ 

1 + P 

Since  and  Qg  are  both  greater  than  1,  (4l)  yields 

2 exp  ' - n log  Qs ) > > expf-  n log  Qg).  (43) 

Thus  the  estimate 


E «=■  v'2 

n 


1 + 6 + */2  6 + 6 ^ ^ 
1 + P 


)■" 


(44) 


is  always  correct  within  a factor  of  4^.  Defining  the  convergence  rate 
to  be 


C = - 6 log  E^/6n,  (45) 

we  can  find  simple  expressions  for  C in  the  limit  6 « 1.  For  p « 6, 
C =»  v'2  6 . When  p » 6,  then  C 6/p.  This  shows  that  convergence 
slows  down  when  the  eigenvalues  of  the  operator  L have  significant 
imaginary  parts. 

Application  to  a Model  Problem 

As  a test  of  the  method  (35)  - (34)  and  the  error  estimate  (37) 
we  shall  solve  numerically  the  two  dimensional  equation 

+ A x,y)*V^  = Six,y),  (46) 

on  a rectangular  grid  with  constant  mesh  spacing  subject  to  doubly 
periodic  boundary  conditions.  Equations  of  this  form  arise  in  plasma 
physics  when  charge  neutrality  is  imposed  upon  a plasma  with  a non- 
isotropic electrical  conductivity,  with  4>  being  the  electrostatic 
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potential.  For  an  application  to  the  earth's  ionosphere,  see  McDonald, 
et  al.^^*^  Another  physical  situation  in  which  an  equation  of  the  form 
(h^S)  arises  is  that  of  a steady  state  distribution  of  a tracer  element 
or  pollutant  ^ in  a given  wind  field  subject  to  known  diffusion  and 
production- loss  rates. 

The  mesh  is  assumed  to  be  by  interior  mesh  points.  We  choose 
the  X direction  so  that 

K ^ K . (i+7) 

X y 

We  also  assume  the  mesh  intervals  to  be 


6x  = 6y  = constant  ('*+8) 

We  use  redundant  guard  cells  around  the  perimeter  of  the  mesh  for 
efficient  handling  of  the  boundary  conditions.  This  results  in  a 
complete  mesh  of  (K^  + 2)  x (K^  + 2)  points.  We  use  second  order  5 
point  representation  of  the  derivatives  in  (46): 


^^IJ  ^^I+1,J  '''  ^I-1,J  ^I,J  + 1 


(49) 


+ AX 
+ AY 


- 4 ^j.j)/6x^ 

- ^ 

I+1,J  I-1,J 

I,J+1  “ ^I,J-1 


)/2  6x 
)/2  6x. 


Here  1 and  J are  the  x and  y mesh  point  indices,  and  AX  and  AY  are  the 
X and  y components  of  respectively. 

In  order  to  estimate  the  eigenvalues  of  L it  is  necessary  to 
consider  ^ constant.  Then  the  eigenfunctions  are  complex  exponentials. 
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i = exp 


where  k = 0,  1,  2, — K - 1,  etc.  for  k . 
X ’ ' ' X ’ y 

The  eigenvalues  are 


COS  + cos  ^ 

I K K 

\ X y 


(AX 
6x 


2TTk  AY  2TTk 

3tn  + _ sin  — 

^ y 


Construction  of  the  Ellipse 

In  order  to  show  that  these  eigenvalues  lie  within  an  ellipse  in  the 
complex  plane,  let  us  adopt  the  following  abbreviated  notation: 


c = cos  — — 

X K 


s = sin  — - 

X K 


a = AX6x,  etc.  for  y. 

X ’ ' 


Let  ? = c + c 

X y 


T]'  = a s + a s . 
XX  y y 


Let  us  consider  the  function 


'.V)  - - (-f-f. 


where  a = (a^  + 

1 X y 
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and  9 is  a tuning  parameter  to  be  determined  later.  Let  us  choose 

a and  a so  as  to  maximize  (55)-  Then 
X y 

h ^ 1/4  (c  + c + r0 f + s^)  (57) 

'xy  xy 

= 1 - 1/4  (c  - C )^  - (1/2  - r9)(s^  + s^), 


where  r = 1 . (58) 

We  shall  carry  r as  a variable  in  the  analysis,  since,  if  we  had 
averaged  h over  all  possible  orientations  of  the  vector  we  would 
have  had  r = 1/2.  It  may  not  be  necessary  in  every  case  to  constrain 
the  method  to  handle  a few  "unfortunate"  orientations  of  A (which 
happen  to  be  resonant  directions  for  the  lowest  few  modes  in  the  system) 
In  general,  A will  vary  in  space,  so  any  destabilizing  resonance  will 
be  localized. 

Let  us  remove  from  consideration  the  modes  (k^,k^)  = (0,0)  and 

K^/2,  Ky/2).  These  are  the  mean  value  and  odd-even  mode  components 

of  The  mean  value  has  no  effect  upon  (47),  and  the  odd-even  mode 

can  be  removed  straightforwardly.  Then  the  lowest  mode  of  consequence 

is  (k  ,k  ) = (1,0),  for  which 
x'  y ' > ' ^ 


s 2TT/K  , K » 1, 

X x'  X ’ 


,59) 


s = 0 . 

y 

To  lowest  order  in  s^,  (57)  gives 

h = 1 - a^(l/2  - r9  ),  where  (60 ) 

max 
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(61) 


We  shall  retain  a as  a variable  since  it  contains  all  relevant  boundary 
condition  information.  Let  us  now  express  S'*  and  Tj'  in  terms  of  the 
real  and  Imaginary  parts  of  in  <52).  Then  ('55)  becomes 


^ h 

max 


(62) 


Comparing  (62)  and  (15),  we  find  the  ellipse  parameters 

a = 4 h^/^  /6x^  (65^ 

max 

b = 4/6 x^ 

c = 9-1/2  -^1/2  /6x2. 

max 

Note  we  have  removed  inconsequential  minus  signs  from  (63).  From 

(39)  and  (63) 


6 = h"i/2  - 1 (64^ 

max 

« 1/4  o2(  1 . 2 re  ), 


Also 

c/a  - P - 1/4  e"i/2. 


(65) 
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Maximizing  the  Convergence  Rate 


The  expressions  (6h)  and  (65)  for  6 and  P allow  the  error  estimate 
(44)  to  be  expressed  in  terms  of  6,  Maximizing  the  convergence  rate 
then  is  equivalent  to  maximizing  Q2  in  (42)  with  respect  to  9 . A 
straightforward  approach  of  setting  BQg/hQ  = 0 requires  solution  of  a 
seventh  degree  polynomial,  A much  more  satisfactory  approach  emerges 
when  one  replaces  (42)  with  an  approximate  expression  valid  for 
5 « 1 ; 

Qo  «=■  1 + v'2  6 + p‘’  - p . (66  ) 

Equation  (66)  gives  correct  results  to  lowest  order  in  6 for  all 
P ^ 0.  Denoting  differentiation  with  respect  to  6 by  prime,  the 
condition  for  maximum  in  (66)  is 


2 6p  '^  = 6 + 2 P6  'p 

Substituting  from  (64)  and  (65),  we  have 


(67) 


r''9^  + n^(5  r9  - 1/2)  = 0, 

where 


(68) 


M- 


a ^1/2 
4 a 

I a1  K 6x 

' X 

8tt 


ri/2 


(69) 


Equations  169)  reflect  the  general  result  that  p and  thus  6 are 
independent  of  the  number  of  cells  into  which  the  rectangular  region 
is  partitioned,  provided  6x  = 6 y,  and  K^,  » 2'^. 
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We  know  that  (68)  has  only  one  real  root,  since  the  derivative 
of  (681  is  never  zero  for  real  0. 

This  root  is 

r0  = I (v'l6  Vi""  + 1 + - (-716  u.'"  + 1 - iP/^  . 

The  solution  to  (68)  has  the  limiting  forms 


p.  « 1 


(71) 


One  can  show  fairly  easily  that  6 varies  monotonically  with  p between 
these  two  limits.  This  behavior  is  verified  in  Fig.  5,  in  which  r6 
is  plotted  vs  M'. 

Convergence  Rates  for  Small  and  Large  | A 1 

Taking  into  consideration  (42),  (44),  (45),  and  66),  we  have  the 
approximate  convergence  rate 


c « ■j2~rTT'-  - 0.  (72) 

We  can  use  (64),  (65  6 nnd  (70 ) to  express  {'J2^  in  terms  of  p..  The 
result  is 

C ctA/2  C (p.),  where  (73) 

o 

C^(p)  = |vi/2  - r0  + p.‘7r0  - p,//r^  j . (7^) 


19 


Recall  that  a = 2tt/K^  for  doubly  periodic  boundaries.  The  dependence 
of  upon  k is  shown  in  Fig.  Also  plotted  in  Fig.  A (dashed  line) 


is  the  approximation 


C^(h.)  (9.5  M.  + 1)’^ 


Using  (71  )>  we  find  the  following  limiting  forms  for  C: 


(75) 


C <=«  c/^2  (1-3  ' u.^/2)i/^j,  p.  « 1 

C ^ ap  ^ /oTo  , P ^ 1 . 


'76) 


This  shows  that  as  p increases  from  0 to  1,  the  convergence  rate  drops 
by  approximately  a factor  of  10. 

Removal  of  Long  Vlavelength  Error  on  Coarse  Grid 

We  can  see  from  (61)  and  (75  that  the  number  of  iterations 
required  to  reach  a certain  level  of  error  reduction  is  proportional 
to  K^.  Thus  the  total  number  of  operations  required  is  proportional 
to  Therefore  it  is  advantageous  to  correct  the  long  wavelength 

components  of  the  solution  on  a coarse  mesh,  and  interpolate  onto  a 
fine  mesh  before  completing  the  solution.  The  Interpolation  introduces 
truncation  error  into  the  solution,  so  that  the  error  estimates  derived 
earlier  may  not  hold  on  the  large  mesh.  For  this  reason,  it  is  best 
to  alternate  between  large  and  small  meshes  more  than  once,  attempting 
only  a modest  error  reductions  with  each  pass.  In  practice,  use  of  a 
coarse  grid  becomes  important  for  > 6k.  For  a grid  of  128  x 128 
interior  points,  a coarse  grid  of  32  x 32  Interior  points  is  used. 
First,  the  residual  ’ - S is  extracted  on  the  fine  grid.  Then  this 


residual  and  A are  defined  on  the  coarse  grid  using  block  averages  of 
l6  fine  grid  points  per  coarse  grid  point.  The  coarse  grid  potential 
correction  is  initialized  to  zero  and  a large  number  of  Iterations  are 
performed.  Typically  we  use  three  times  as  many  iterations  on  the 
coarse  grid  as  on  the  fine  grid.  Then  the  potential  correction  is 
defined  on  the  fine  grid  using  bilinear  interpolation,  and  the  result 
is  subtracted  from  Then  iterations  are  performed  on  the  fine 

grid,  making  no  attempt  to  make  significant  improvements  in  the  lowest 
5 to  10  modes  of  the  solution.  That  is,  we  arbitrarily  increase  o 
in  61^  by  a factor  of  5 to  10.  This  results  in  Improved  convergence 
of  the  higher  modes.  The  effectiveness  of  this  fine  grid-coarse  grid 
approach  may  be  improved  significantly  in  a time  dependent  problem  by 
two-level  extrapolation  for  the  trial  solution  <h^ . 

Numerical  Results 

As  a test  of  the  iterative  procedure  (35)  - (3^)  and  the  convergence 
estimate  (73)  have  solved  (1*6)  on  meshes  of  52  x 32  and  US  x l*-8 
interior  grid  points  without  regrldding;  and  on  a 128  x 128  mesh  with 
a reduced  mesh  of  32  x 32  interior  points.  All  tests  used  doubly 
periodic  boundary  conditions.  The  numerical  convergence  rate  comparisons 
were  carried  out  as  follows.  Arbitrary  forms  were  adopted  for  A and  a 
reference  solution  $ . Then  a source  term  was  generated  from  numeri- 
cally using  the  difference  operation  (49^.  This  source  term  was  used 
in  the  iteration  (53)  - with  the  approximate  solution  being 

Initialized  to  zero.  After  a large  number  N of  Iterations  (usually 
N = 1*0),  a relative  error  E was  defined  to  be  the  root  mean  square 
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residual  of  ,1+6)  divided  by  the  root  mean  square  of  the  source,  S. 

The  average  convergence  rate  was  then  taken  to  be  - log  E/N. 

In  all  cases  convergence  was  fast  enough  to  be  consistent  with 
73  . and  in  most  cases  was  faster  than  (73)*  Th®  one  case  in  which 
convergence  was  just  equal  to  75)  fot  sH  h-  was  for  ~ sin 
c!  n i/k  , and  A = constant,  A =0.  This  is  consistent  with  the 

XX  y 

discussion  following  (58).  For  determining  - in  '69),  we  used  the 
root  mean  square  value  for  , and  r = 1. 

Figure  5 shows  convergence  rates  per  second  of  computer  time, 

- log  E/t,  obtained  using  a two  pipe  Texas  Instruments  ASC  to  solve 
the  same  set  of  problems  with  the  Chebychev  explicit  method  CE ) 
upper  curve  and  with  an  alternating  direction  implicit  method  ADI ) 
lower  curve  . Thcoe  comparisons  were  made  on  a 5O  x 50  mesh  using 
data  from  the  EJET  plasma  turbulence  code/’^'  The  ADI  solution  used 
logarithmically  spaced  iteration  parameters'^  and  a partially  vectorized 
tridiagonal  solver.'  Although  ADI  converged  faster  per  iteration,  the 
limited  vectorizabi llty  of  the  method  resulted  in  slow  execution.  Thus 
the  result  in  Fig.  F is  machine  dependent.  The  execution  times  per 
iteration  were  25-8  milliseconds  for  ADI  and  1.32  milliseconds  for  t'>e 
Chebychev  explicit  method.  To  illustrate  the  computational  efficiency 
of  the  explicit  method,  the  lower  limit  on  execution  time  per  iteration 
neglecting  overhead  and  boundary  value  resetting)  would  be  (4H  x W 
Interior  points)  x (11  operations  per  point)  x (40  nanoseconds  per 
operation)  = 1.01  milliseconds.  Thus  the  explicit  method  runs  at  Jo'i 
efficiency  on  a modest  50  50  mesh.  ADI  could  be  made  more  competitive 
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by  the  use  of  cyclic  reduction  rather  than  tridiagonal  solution  in  the 
integration  direction,  but  it  is  unlikely  that  any  improvement  would 
bring  the  convergence  rate  per  second  up  to  that  of  the  explicit  method. 

In  addition,  the  boundary  conditions  are  much  easier  to  change  in  the 
explicit  method  than  in  the  implicit  one. 
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Fig.  2 — Krror  estimate  E,,  versus  n. 
Curves  are  labelled  with  log]Q  R. 
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Fig.  4 — Normalized  convergence  rate  (74)  (solid); 
approximation  (75)  (dashed) 


28 


I 


Fig.  5 — Computational  efficiencies  of  the  Chehychev  explicit  method 
(CK)  and  alternating  direction  implicit  method  (ADI)  when  applied  to 
a senes  of  numerical  tests  in  which  lA  I and  thus  p vary  from  case  to 
case.  Plotted  is  the  logarithm  of  the  error  reduction  factor  after  many 
Iterations,  divided  hy  the  elapsed  computer  time  in  seconds. 


I 


i 


29 


DISTRIIUTION  LIST 

DIRECTOR 

Defense  Advanced  Rsch  Proj  Agency 
Architect  Building 
ll*00  Wilson  Blvd. 

Arlington,  VA  22209 

ATTN:  Strategic  Tech  Office 

Defense  Communication  Engineer  Center 
1660  Wiehle  Avenue 
Res  ton,  VA  22090 

ATTN:  CODE  R410  W.  D.  Dehart 

DIRECTOR 

Defense  Communications  Agency 
Washington,  D.  C.  20305 
ATTN : CODE  96O 

ATTN : CODE  480 

Defense  Documentation  Center 
Cameron  Station 

Alexandria,  VA  22514  12  copies  (if  open  publication) 

ATTN:  TC  2 copies  (if  otherwise) 

DIRECTOR 

Defense  Intelligence  Agency 
Washington,  D.  C.  20501 

ATTN:  W.  Wittig  DC-7D 

ATTN:  DT-lB 


DIRECTOR 

Defense  Nuclear  Agency 

Washington, 

D.  C.  20305 

ATTN: 

STS  I Archives 

ATTN: 

STVL 

ATTN: 

STTL  Tech  Library 

2 copies 

ATTN: 

DOST 

ATTN: 

RAAE 

3 copies 

DIR  OF  DEFENSE  RSCH  & ENGINEERING 
Washington,  D.  C.  20501 

ATTN:  DD/S&SS  John  B.  Walsh 

ATTN:  OAD/EPS  LTC  W.  A.  Whitaker 


1 


COMMANDER 
Field  Command 
Defense  Nuclear  Agency 
Kir t land  AFB,  NM  87115 
ATTN ; FCPR 

ATTN:  FCPR  COL.  John  P.  Hill 

Interservice  Nuclear  Weapons  School 
Kirtland  AFB,  NM  87II5 

ATTN:  Document  Control 

DIRECTOR 

Joint  Strat  TGT  Planning  Staff  JCs 
Offutt  AFB 
Omaha,  NB  68113 
ATTN:  JLTW-2 

ATTN:  JPST  G.  D.  Burton 

ATTN:  >TPST  MAJ.  J.  S.  Green 

CHIEF 

Livermore  Division  FLD  Command  DNA 
Lawrence  Livermore  Laboratory 
P.  0.  Box  808 
Livermore,  CA  9^550 
ATTN:  FCPRL 

COMMANDER 

National  Military  Comd  Sys  Support  Ctr 
Pentagon 

Washington,  D.  C.  20501 
ATTN:  8211 

ATTN:  DP  Director  for  CSPO 

OJCS/J-3 
The  Pentagon 

Washington,  D.  C.  20501 

ATTN:  J-5  OPS  ANAL  BR.  COL.  Longsberry 

ATTN:  J-6 

DIRECTOR 

Telecommunications  & Comd  & Con  Sys 
Washington,  D.  C.  20301 

ATTN:  ASST  DIR  Info  & Space  Sys. 

ATTN:  DEP  ASST.  SEC  Sys. 

DIRECTOR 

National  Security  Agency 
Ft.  George  G.  Meade,  Md.  20755 
ATTN:  Wl4  Pat  Clark 

ATTN : Frank  Leonard 


2 


Weapons  Systems  Evaluation  Group 
^00  Army  Navy  Drive 
Arlington,  VA  22202 

ATTN:  D0CU>ENT  CONTROL 

DIRECTOR 

BMD  Advanced  Tech  Ctr. 

Huntsville  Office 
P.  0.  Box  1500 
Huntsville,  AL  55807 

ATTN:  ATC-T  Melvin  T.  Capps 

MANAGER 

Bl©  Program  Office 
1500  Wilson  Blvd. 

Arlington,  VA  22209 

ATTN:  Plans  Division 

ATTN:  DACS-BMS  Julian  Davidson 

COMMANDER 

Harry  Diamond  Laboratories 
2800  Powser  Mill  Road 
Adelphi,  Md.  20783 
ATTN : AMXDO-NP 

COMMANDER 

TRASANA 

White  Sands  Mlssle  Range,  NM  88002 
ATTN:  EAB 

DIRECTOR 

U.S.  Army  Ballistic  Research  Labs 
Aberdeen  Proving  Ground,  Md.  21003 
ATTN:  AM-CA  Franklin  E.  Niles 

U.S,  Army  Communications  CMD 
C-E  Services  Division 
Pentagon  Rm.  2d513 
Washington,  D.  C.  20310 
ATTN:  CEAD 


COMMANDER 

U.S.  Army  Material  Command 
5001  Eisenhower  Avenue 
Alexandria,  VA  22533 

ATTN:  AMCRD-WN-RE  John  F.  Corrigan 

Commander 

U.S.  Army  Electronics  Command 
Fort  Monmouth,  N.  J.  07703 
ATTN:  AMSEL-TL-ENV  Hans  A.  Bomke 


5 


COMMANDER 

U.S.  Army  Material  Command 
Foreign  and  Scientific  Tech  Center 
290  7th  St.  N.  E. 

Charlottesville,  VA  22901 
ATTN:  r.  A.  Crowley 

ATTN:  R.  Jones 

COMMANDER 

U.S.  Army  Missile  Command 
Redstone  Arsenal 
Huntsville,  AL  35809 

ATTN:  ASMI-YTT  W.  G.  Preussel,  Jr. 

COMMANDER 

U.S.  Army  Nuclear  Agency 
Fort  Bliss,  TX  79916 
ATTN:  CD INS -E 

ATTN:  USANUA-W.  J.  Herbert 

CHIEF  OF  NAVAL  RESEARCH 
Department  of  the  Navy 
Arlington,  VA  22217 

ATTN:  CODE  h6k  Jacob  L.  Warner 

ATTN : CODE  h6k 

COMMANDER 

Naval  Air  Systems  Command 
Headquarters 

Washington,  D.  C.  21550 
ATTN:  AIR  558I 

COMMANDER 

Naval  Electronics  Systems  Command 
Naval  Electronics  Systems  CMD  HQS 
Washington,  D.  C.  2056O 

ATTN:  NAVALEX  034  T.  Barry  Hughes 

ATTN:  PME  I06-I  Satellite  Conin  Project  Off 

ATTN:  John  E.  Doncarlos 

ATTN:  PME  II7 

COMMANDER 

Naval  Electronics  Laboratory  Center 
San  Diego,  CA  92152 

ATTN:  William  F.  Moler 

ATTN:  CODE  2200  1 Verne  E.  Hildebrand 

ATTN:  R.  Eastman 


4 


COMMANDING  OFFICER 
Naval  Intelligence  Support  CTR 
4301  Suit  land  Road,  Bldg.  5 
Washington,  D.  C.  20390 

ATTN:  Mr.  Dubbin  Stic  12 

DIRECTOR 

Naval  Research  Laboratory 
Washington,  D.  C.  20375 

ATTN:  HDQ  COMM  DIR  Bruce  Wald 

ATTN:  CODE  ^h60  Radio  Propagation  BR 

ATTN:  CODE  7127,  Charles  Y. Johnson 

ATTN:  CODE  7701,  Jack  D.  Brown 

ATTN:  CODE  7700,  Division  Superintendent,  25  copies  (if  open  publication) 

1 copy  (if  otherwise) 

ATTN:  CODE  7750,  Branch  Head  I50  copies  (if  open  publication) 

1 copy  (if  otherwise) 


COMMANDING  OFFICER 

Naval  Space  Surveillance  System 

Dahlgren,  VA  22448 

ATTN:  Capt.  J.  H.  Burton 

COMMANDER 

Naval  Surface  Weapons  Center 
White  Oak,  Silver  Spring,  Md.  209IO 
ATTN:  CODE  730  Tech.  Lib. 

ATTN:  CODE  1224  Navy  Nuc  Prgms  Off 

DIRECTOR 

Strategic  Systems  Project  Office 
Navy  Department 
Washington,  D.  C.  2037^ 

ATTN:  NSP-2141 


COMMANDER 

ADC/AD 

ENT  AFB  CO  50912 
ATTi; : ADDA 

AF  Cambridge  Rsch  Labs,  AFSC 
L.  G.  Hanscom  Field 
Bedford,  MA  01750 

ATTN:  LKB  Kenneth  S.  W.  Champion 

ATTN:  OPR  Harvey  P.  Gauvin 

ATTN:  OPR  James  C.  Ulwick 

AF  Weapons  Laboratory,  AFSC 
Kirtland  AF’B,  NM  87117 


ATTN: 

DYT  LT  Mark  A.  Fry 

ATTN: 

DYT  CAPT  Wittwer 

ATTN: 

John  M.  Kamm  SAS 

ATTN: 

SUL 

5 


Ai'TAC 

Patrick  Ai’B, 
ATTN: 
ATTN; 
ATTN: 


FL  3292^ 

TF  MAJ.  E.  Hines 
TF/CAPT.  Wiley 
TN 


Air  Force  Avionics  Laboratory,  AFSC 
Wright-Patterson  AFB,  OH  U5U33 

ATTN:  AFAL  AVWE  Wade  T.  Hunt 

Assistant  Chief  of  Staff 
Studies  and  Analysis 
Headquarters,  U.  S.  Air  Force 
Washington,  D.  C.  20330 

Headquarters 

Electronic  Systems  Division,  (AFSC) 
L.  G.  Hanscom  Field 
Bedford,  MA  01730 

ATTIi:  XFF:  LT.  Michaels 

ATTN:  LTC  J.  Morin  CDEF  XRC 

ATTN:  YSEV 

COMMANDER 

Foreign  Technology  Division,  AFSC 
Wright-Patterson  AFB,  OH  U5U33 
ATTN:  TD-BTA,  LTBRARY 

HQ,  USAF/RD 

V/ashington,  D.  C.  20330 
ATTIi : RDQ 

COICIANDER 

Rome  Air  DeveT''pment  Center,  fC 
Griffiss,  AFB,  N.  U.  13l+i+0 
ATTN:  EMTIiD  Doc  Library 


Ca-nVITiER  nc  CHIEF 
Strategic  Air  Command 
Offutt  /F'B,  NB  68113 

ATTN:  XPFS  fiAJ.  Brian  G.  Stephan 

544  lES 

Offutt  AFB,  NB  68113 

ATTN:  RDPO  LT.  Alan  B.  Merrill 
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Los  Alamos  Scientific  Laboratory 

P.  0.  Box  1663 

Los  Alamos,  in-1  8?5^^ 

ATTN:  DOC  CON  for  Eric  Lindman 

ATTN:  DOC  CON  for  H.  M.  Peek 

ATTN:  DOC  CON  for  R . F.  Taschek 

Sandia  Laboratories 
P.  0.  Box  5800 
Albuquerque,  NM  87115 

ATTN:  DOC  CON  for  A.  Dean  Thronbrough 

ATTN:  DOC  CON  for  W.  D.  Brown 

ATTIi:  DOC  C^N  for  D.  A.  Dahlgren,  ORG  1722 

ATTN:  DOC  CON  for  J.  P.  Martin,  ORG  1752 

University  of  California 
Lawrence  Livermore  Laboratory 
P.  0.  Box  808 
Livermore,  CA  9^550 

ATTN:  Tech  Info  Dept  L-3 

Department  of  Commerce 
National  Oceanic  & Atmospheric  Admin. 
Environmental  Research  Laboratories 
Bo’ulder,  CO  80302 

ATTN:  Jospeh  H.  Pope 

ATTN:  C.  L.  Rufenach 

Departeent  of  Commerce 
Office  of  Telecomm.iunications 
Institute  for  Telecom  Science 
Boulder,  CC  80302 

ATTN:  Glenn  Falcon 

ATT?:  V.’illiam  F.  Utlaut 

A iTI. : G , Reed 

ATT??:  L.  Berry 

Department  of  Transportation 
Transportation  Rsch.  System  Center 
Kendall  Square 
Cambridge,  M/\.  02142 

ATTN:  TER  G.  Harowles 

NASA 

Goddard  Space  Flight  Center 
Greenbelt,  iD  20771 

ATTN:  CODE  750  T.  Golden 
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NASA 

600  Independence  Ave . , S.  W. 
Washington,  D.  C.  205^ 

ATTN ; H . Dub  in 

Aerospace  Corporation 
P.  0.  Box  92957 


Los  Angeles,  CA  90009 


ATTN: 

T.  M.  Salmi 

ATTN; 

S . P . Bower 

ATTN; 

V . Josephson 

ATTN: 

SMFA  for  PWW 

Analytical  Systems  Engineering  Corporation 
S Old  Concord  Rd. 

Burlington,  MA  O8IO5 
ATTN:  Radio  Sciences 

Avco-Everett  Research  Laboratory,  Inc. 

2585  Reverse  Beach  Parkway 
Everett,  MA  021^9 

ATTN:  Richard  M.  Patrick 

Boeing  Company,  The 
P.  0.  Box  5707 
Seattle,  WA  98124 
ATTN : D . Murray 

ATTN:  Glen  Keister 

Brown  Engineering  Company,  Inc. 

Cummings  Research  Park 
Huntsville,  AL  55807 

ATTN:  David  Lambert  MS  I8 

California,  University  of  - Sen 
Sfccnrit^,  Cffice  --iv 
La  Jolla , .’A  ,.'2098 

ATTN:  Henry  G.  Booker 

Caispan  CoT-noration 
P.  0.  Box  . 5 
Buffalo,  N.  14221 

ATTN:  Romeo  A.  De liber is 

Computer  Sciences  Corporation 
P . 0 . Box  550 
6565  Arlington,  VA  22046 
ATTN:  H.  Blank 
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Comsat  Laboratories 
P.  0.  Box  115 
Clarksburg,  Md . 20734 
ATTN:  R.  R.  Taur 

Cornell  University 

Department  of  Electrical  Engineering 
Ithaca,  N.  Y.  14850 

ATTN:  D.  T.  Farley,  Jr. 

ESL,  INC. 

495  Java  Drive 
Sunnyvale,  CA  94o86 

ATTN:  J,  Roberts 

ATTN:  V.  L. Mower 

ATTN : James  Marshall 

ATTN:  R.  K.  Stevens 

General  Electric  Company 
Ten^jo-Center  for  Advanced  Studies 
816  State  Street 
Santa  Barbara,  CA  93102 
ATTN:  Don  Chandler 

ATTN:  DASIAC 

ATTN:  Tim  Stephens 

General  Electric  Company 
P.  0.  Box  1122 
Syracijse,  N.  Y.  13201 

ATTN:  F.  A.  Reibert 

General  Research  Corporation 
P.  0.  Box  3587 
Santa  Barbara,  CA  93105 
ATTN:  John  Ise,  Jr. 

Gecpy-.: ical  Institute 
University  01  /JLaska 
F-ircanks,  AK  99701 

AIT.'.:  Technical  ibrary 

ATir:  Neil  Bro'.'m 

AtTN:  T.  N.  Davis 

GTE  Syl '.a,  Inc. 

189  B Street 

Needham  Heights,  M\  02194 
aTTN:  Miars^iUl  Cross 

HRB-Singer,  Inc. 

Science  Park,  Science  Park  Road 

P.  0.  Box  60 

State  College,  PA  I68OI 

ATTN:  larry  Feathers 
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UUnois, 

133  Street 

Bo7  Sovtb  .^B20 

CY^.ar^paipn,  ^.^aiy^es 

! for  Defense  Knni.^ 

instit^i^®  Dri'^ 

400  S 22202 

BrlinS^^’  Ernest  Bauer 
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. Teiegrant^  Corp 
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Marietta  Corporation 
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Denver  D 
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'>*'"'=VS:  special  Proj:: 


Maxwell  Tjaborator ies , Inc. 

92U4  Balboa  Avenue 
San  Diego,  CA  ^*2123 

ATTI.':  A.  J.  Stormon 

ATTIi:  V.  Fargo 

ATTII;  A.  N.  Rostocker 

McDonnell  Douglas  Corporation 
5301  Bolsa  Avenue 
Huntington  Beach,  CA  92647 
ATTN : J.  Moule 

ATTN:  N.  Harris 

ATTN:  Tech  Library  Services 

Mission  Research  Corporation 
735  State  Street 
Santa  Barbara,  CA  93101 
ATTN:  Ralph  Kilb 

ATTN:  R.  Hendrick 

ATTN:  Conrad  L.  Longmire 

ATTN:  R.  E.  Rosenthal 

ATTN:  R.  Bogusch 

ATTI'I:  David  Sowle 

ATTN:  M.  Scheibe 

ATTii;  P.  Fischer 

Mitre  Corporation,  The 
Route  52  and  Middlesex  Turnpike 
P.  0.  Box  208 
Bedford,  Ilk  01730 

ATT.::  S.  A.  Morin  M/3 

A'TT;:  Chief  Scientist  Sen 

-ATT.':  G.  :iarding 

Mitre  Corporation,  The 
'..estgate  Research  lark 
1820  Dolley  Madison  Blvd. 

McLean,  VA  22101 

ATT;:  ll-.-n  Scimeider 

North  C-jrolina  State  Univ.  at  Raleigh 
North  Carolina  State  Univ.  Carpus 
Raleigh,  N.  c.  27606 

ATTN:  SEC  OFFICER  for  Walter  A.  Flood 

Pacific -Sierra  Research  Corp. 

1456  Cloverfield  Blvd. 

Santa  Monica,  CA  9-4-4 

ATTN:  E.  C.  Field,  Jr. 
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Philco-Ford  Corporation 
'.-.’estern  Developrr.ent  Laboratories  Div. 
3939  Fabian  '•'ay 
Falo  Alto,  CA  9^*303 

ATT:I:  J.  T.  Mattingley  MS  X22 

Photoretrics , Inc. 

UU3  Marrett  Road 
Lexington,  MA  02173 

ATTN;  Irving  L.  Kofsky 

Physical  Dynamics,  Inc. 

P.  0.  Box  1069 
Berkeley,  CA  9^701 

ATTT^;  Joseph  B.  V/orkman 

R 8c  D Associates 
P.  0.  Box  3580 
Santa  Monica,  CA  90^02 


ATTN; 

Forest  Gilmore 

ATTN; 

Richard  Latter 

ATTI.’; 

Robert  E.  Lelevier 

ATTN; 

V/illiam  B.  Wright, 

Band  Corporation,  The 

1700  Mai 

n Street 

Santa  Monica,  CA  90^C6 

A X Xl'l  .* 

Cullen  Crain 

Science 

Applications,  Inc. 

?.  C.  rc 

;^*2351 

La  Leila 

, CA  92038 

J * 

Eaniel  A Hamlin 

- - . ♦ t 

1 • s 

Z.  Staraker 

Stanford 

Research  Institute 

333  Ravenswood  Avenue 

■■•=nlo  Pa 

rk,  CA  '>1025 

a:  IN; 

L.  L.  Cobb 

attt:; 

Walter  G.  Chestnut 

ATTN; 

Donald  Neilson 

ATTN; 

David  A.  Johnson 

ATTN; 

Charles  L.  Rino 

att:;-. 

E.  J.  Fremouw 

Air:; 

Dale  H.  Davis 

TR’i  Systems  Group 
One  Space  Park 
Pedondo  Beach,  CA 
ATTN;  J.  V/.  Lowry 
ATTN;  P.  H.  Katsos 
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Advanced  Research  Projects  Agency  (ARPA) 
Strategic  Technology  Office 
Arlington,  Virginia 

Capt.  Donald  M.  LeVine 

Naval  Research  Laboratory 
V;ashington,  D.  C.  20575 


Dr.  P.  Mange 
Dr.  E.  Peterkin 
Dr.  R.  Meier 
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Dr.  J.  Goodman  --  Code  7950 

Science  .applications,  Inc. 
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F.  Donnelly 
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Bedford,  Mass.  0l7p0 
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Dr.  W.  Swider 
Mrs.  R.  Srgal>Ti 
Dr.  J.  M.  Forbes 
Dr.  T.  .1.  Keneshoa 

Office  of  Naval  Research 
80Z  North  Quincy  Street 
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Dr.  J.  G.  Dardis 
Dr.  H-  Mullaney 
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Dr. 
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U.  S.  Army  Aberdeen  Research  and  Development  Center 
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Naval  A.ir  Systems  Command 
Department  of  the  Navy 
Washington,  D.  C.  20^G") 

Dr.  T.  Czuba 


Harvari  University 
Harvard  Square 
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Dr.  F.  V.  Coroniti 
Dr.  C.  Kennel 


14 


University  of  California,  San  Diego 
P.  0.  Box  109 
La  Jolla,  California 


Dr.  P. 

M.  Banks 

Utah 

State 

University 

4 th 

N . and 

8th  E.  Strei 

Loga 

n,  Utah 

84522 

Dr.  R. 

Harris 

Dr.  V. 

Peterson 

Dr . R . 

Megill 

Cornell  University 
Ithaca,  New  York  14850 
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